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Abstract 

If an evaporating black hole does not settle down to a non radiating remnant, a 
description by a semi classical Lorentz metric must contain either a naked singularity or 
what we call a thunderbolt, a singularity that spreads out to infinity on a spacelike or null 
path. We investigate this question in the context of various two dimensional models that 
have been proposed. We find that if the semi classical equations have an extra symmetry 
that make them solvable in closed form, they seem to predict naked singularities but 
numerical calculations indicate that more general semi classical equations, such as the 
original CGHS ones give rise to thunderbolts. We therefore expect that the semi classical 
approximation in four dimensions will lead to thunderbolts. We interpret the prediction of 
thunderbolts as indicating that the semi classical approximation breaks down at the end 
point of black hole evaporation, and we would expect that a full quantum treatment would 
replace the thunderbolt with a burst of high energy particles. The energy in such a burst 
would be too small to account for the observed gamma ray bursts. 
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1 Introduction 

It has been known for some time that classical general relativity predicts singularities 
in gravitational collapse. At the singularities, the Einstein equations will not be defined. 
Thus there will be a limit as to how far into the future one can predict spacetime. However, 
it seems that singularities formed in gravitational collapse always occur in regions that are 
hidden from infinity by an event horizon, so the breakdown of the Einstein equations at 
the singularity does not affect our ability to predict the future in the asymptotic region of 
space. This assumption that the singularities are hidden is known as the Cosmic Censorship 
Hypothesis and is fundamental to all the work that has been done on black holes. It 
remains unproven but it is almost certainly true for classical general relativity with a 
suitable definition of a singularity that is so bad it can't be smoothed out or continued 
through. 

On the other hand, in the semi classical approximation to quantum gravity a black 
hole formed in a gravitational collapse will emit thermal radiation and evaporate slowly. 
If the black hole has a charge that is coupled to a long range field and which can't be 
radiated, such as a magnetic charge, it may be able to settle down to a non radiating 
state such as the extreme Reissner-Nordstr0m solution. But for black holes without such 
a charge, there are no zero temperature classical solutions they can settle down to. One 
might suppose they settled down to some stable or semi stable remnant that was not a 
classical solution but was maintained by quantum effects. However, quite apart from the 
fact that there is nothing very obvious to stabilize such remnants, their existence would 
create severe problems. If they had a mass of the order of the Planck mass, one might have 
expected that there would be more than the cosmological critical density of the remains 
of black holes formed in the very early universe. While if they had zero mass, they would 
lead to infinite degeneracy of the vacuum state. 

The most natural assumption would seem to be that black holes without a conserved 
charge disappear completely. To suppose that black holes could be formed but never 
disappear would violate CPT unless there were also a separate species of white holes 
which would have existed from the beginning of the universe. On the other hand, if black 
holes disappear completely, black and white holes can be different aspects of the same 
objects, which would be an aesthetically satisfying solution to the CPT problem. Holes 
would be called black when they were large and classical, and not radiating much, but 
they would be called white when the quantum emission was the dominant process. 

If black holes disappear completely, this can not be described by a Lorentzian metric 
without some sort of naked singularity, or what would be even worse, a region of closed 
time like curves. Spreading out from the naked singularity or region of chronology violation 
would be a Cauchy horizon. Beyond this horizon the semi classical equations would not 
uniquely specify the solution, but one would hope that it would determined by a full 
quantum treatment, though maybe with loss of quantum coherence. Otherwise, we could 
be in for a surprise every time a black hole on our past light cone evaporates. 

Within the context of the semi classical approximation there is however an alternative 
to a naked singularity that has not received much attention. We shall call it a thunderbolt. 
It is a singularity that spreads out to infinity on a space like or null path. It is not a 
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naked singularity because you don't see it coming until it hits you and wipes you out. 
It would mean that the semi classical equations could not only not be evolved uniquely 
(as with a naked singularity), but they could not be evolved at all more than a finite 
distance into the future. If the thunderbolt was null, one could regard it as the singular 
Cauchy horizon produced by some would-be naked singularity. This would be like what 
is believed to happen to the inner Cauchy horizons of classical black holes under generic 
perturbations. One might therefore expect that although the semi classical equations 
could lead to naked singularities in special situations, one would get a thunderbolt if one 
perturbed the equations or the initial data slightly. 

If the semi classical equations were to predict a thunderbolt singularity as the end 
point of black hole evaporation, one would have to conclude that the singularity would be 
softened and smeared out by quantum effects because surely many black holes must have 
evaporated in the past, and yet we survived. Nevertheless, if the semi classical equations 
predict thunderbolts, this might indicate that something fairly dramatic happens in the 
full quantum theory. 

In four dimensions, the one loop corrections are quadratic in the curvature. This 
means that the semi classical equations including one loop back reaction are fourth order 
and have unphysical runaway solutions. It is therefore hard to use them to decide whether 
the evaporation of black holes leads to naked singularities or thunderbolts. On the other 
hand, the the one loop corrections in two dimensions are proportional to the curvature 
scalar. This means that the semi classical equations are second order even when the back 
reaction is taken into account. It should therefore be possible to decide what they predict 
as the outcome of black hole evaporation. Hopefully, this will give an indication of what 
might happen in four dimensions. 

In two dimensions the Einstein Hilbert Lagrangian R is a divergence. This means 
that to get a non trivial interaction with the metric, one has to multiply the Einstein 
Hilbert term by a function of a dilaton field (p. An interesting model in which the metric is 
coupled to a dilaton field and N minimal scalars has been proposed by Callan, Giddings, 
Harvey and Strominger [1], (henceforth referred to as CGHS). In the classical version of 
this theory one can form a black hole by sending in a wave of one of the scalar fields from 
the asymptotic region. Quantum field theory on this classical black hole background then 
shows that the black hole will radiate thermally in each of the fields. Presumably this 
means that the black holes will evaporate but a full quantum treatment of the problem 
seems too difficult even in this simple theory. However, Callan et al suggested that in 
the large N limit, one could neglect ghosts and quantum fluctuations of the metric and 
dilaton in comparison with those of the scalar fields. The effective action arising from 
the scalar quantum loops would be completely determined by the trace anomaly and the 
conservation equations together with boundary conditions. One could therefore add it 
to the classical action for the metric and dilaton fields and obtain a set of semi classical 
hyperbolic differential equations for the metric and dilaton. 

Even these relatively simple equations have not been solved in closed form. Callan 
et al hoped that the result of including the action of the scalar loops would be to cause a 
black hole to evaporate completely without any singularity and tend at late times to the 
linear dilaton solution, which is the analogue of Minkowski space, and which is the natural 
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candidate for a ground state. However, later work showed that there was necessarily a 
singularity, and that the solution could not settle down to a static state in which the 
singularity remained hidden behind an event horizon. 

These results presumably indicate that the semi classical equations lead either to a 
naked singularity or a thunderbolt. But which? The original semi classical equations 
proposed by CGHS do not seem to admit closed form solutions. Various authors have 
suggests modifications to the semi classical equations that introduce an extra symmetry 
and make the equations solvable in closed form. We shall show the exact solutions have 
naked singularities. However they also continue to emit radiation at a finite rate and 
the mass becomes arbitrarily negative. Such behaviour is presumably unphysical, or at 
least one hopes so. The conservation of energy would lose its practical significance if 
one could have negative mass naked singularities. In one case at least, one could use 
the non uniqueness of the solution after the naked singularity has appeared to cut off 
the analytically continued exact solution at the Cauchy horizon produced by the naked 
singularity and glue on a non radiating solution. This procedure however transforms the 
Cauchy horizon into a thunderbolt singularity, although a fairly mild one. 

In the four dimensional case, the equations don't have symmetries that allow one 
to solve them in closed form. There is thus no reason to expect special properties like 
conformal symmetry in two dimensional models of black holes. We shall therefore investi- 
gate the behaviour of solutions of the original semi classical equations proposed by CGHS 
which we expect to be more typical of the general case. Since these equations do not 
admit solutions in closed form, there seems no alternative but to integrate the equations 
numerically. Fortunately hyperbolic equations in 1+1 dimensions are relatively easy and 
there are reliable and numerically stable routines available. To test their accuracy, we first 
applied them to the equations without back reaction. We obtained excellent agreement 
with the known solution, the Witten two dimensional black hole. Encouraged by this, we 
included the back reaction terms and obtained results that strongly indicate a thunderbolt. 
This supports our view that while naked singularities may occur for certain sets of semi 
classical equations with special symmetries, more general two dimensional models of black 
hole evaporation will exhibit thunderbolts. 

In section 2 the model and the various sets of semi classical equations are described. 
Those with special symmetries that allow exact solutions are shown to lead to naked 
singularities in section 3 while in section 4 the numerical results of integrating more general 
equations are presented. A test is given to distinguish a thunderbolt from an eternal black 
hole. The implications for black holes in four dimensions are discussed in section 5. The 
numerical algorithm used is described in an appendix. 
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2. The semi classical model 



CGHS assume the spacetime contains a dilaton field <fi and N minimally coupled scalar 
fields /j, described by a classical Lagrangian 



2tt 



-g 



-20 



1 N 

(i? + 4(V^) 2 + 4A 2 )--^(V/,) 



i=i 



(2.1) 



where R is the Ricci scalar and A is a coupling constant. 

Any two dimensional spacetime is of course conformally flat, so one can introduce null 
coordinates x ± and write the line element as 



ds 2 = -e 2p dx+ dx 



(2.2) 



CGHS suggested that in the limit of a large number iV of scalar fields fi one could neglect 
the quantum fluctuations of the dilaton and the metric, and treat the back reaction in the 
scalar fields semi classically by adding to the action a trace anomaly term 



—K,d + p d-p. 

CGHS took k = N/12. However taking ghosts into account leads to 

N -24 



K = 



12 



(2.3) 



(2.4) 



in that theory. For consistency with refs [3-5] we henceforth define n by (2.4). Occasionally 
we shall use the earlier value in the form R = N/12; obviously h = k + 2. We shall call the 
theory defined by equations (2.1), (2.3) and (2.4) the original theory. 

Strominger [2] has suggested that the ghosts should be coupled to a different metric. 
This leads to the action of the original theory (with k replaced by k), plus an additional 
term 

2 (d+(/>d-(/> - d+c/>d-p - d+p <9_0 + d+p d_p) . (2.5) 

We shall call this the decoupled ghost theory, though in fact the ghosts are still coupled 
to the geometry, only differently. 

de Alwis [3] and Bilal and Callan [4] have suggested that the cosmological constant 
A 2 term be multiplied by a function D((p) to make the theory conformally invariant where 



D((f>) = ^(l + y) 2 exp 



1 



1 + y 



(2.6) 



and y = y/l — ne 2 ^. We shall call this the conformal theory. It can be solved in closed 
form. 

Another Lagrangian with a special symmetry that has a conserved current j M = 
dfj,(4> — p) has been proposed by Russo, Susskind and Thorlacius [5]. It is the Lagrangian 
of the original theory plus the additional term 

-K(/>d+d-p (2.7) 
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We shall call this the conserved current theory. 

The general solution of the conformal and conserved current theories with an asymp- 
totically flat weak coupling region will be given in section 3. It will be shown they have 
naked singularities for positive k. Here we give the field equations for the two Lagrangians 
without special symmetries, the original and decoupled ghost theories. The evolution 
equations can be written in the form 

d+d-fr =0, (2.8a) 
d+d-p =p- 1 {2d + (t)d-(t) + Y), (2.86) 
d+d-<f) =Qd+d-p, (2.8c) 



where we have introduced the quantities 



in the original theory, and 



P =1 - Ke 2<t> 
Q=l- ifte 2 ^, 

p =1 - he 2 ^ + \ke^ 



(2.9a) 



(2.96) 



Q=l- §«e 2 ^, 
in the decoupled ghost theory. Here 

Y = \^e 2p . (2.10) 
In addition there are two constraint equations. In the original theory they are 

e- 2 * {2d%<j> - 4d +( pd + p) - k [d 2 + p - (d+p) 2 - t + (x+)} =\ J2(d+h) 2 , (2-lla) 

i 

e- 2 * {2d 2 _<P - 4cL0d_p) - k [d 2 _p - (d_p) 2 - t.(x-)] =\ J2(d-fi) 2 , (2.H6) 

i 

where t± are arbitrary functions. They are constraints in the following sense. (2.11a,b) 
need be imposed only on surfaces x~ = const, and x + = const, respectively. They hold 
then throughout the spacetime as a consequence of the evolution equations. The con- 
straints for the decoupled ghost theory involve replacing k by k as well as adding some 
extra terms which vanish when <p = p. Since we only impose the constraints on the initial 
surfaces where we may also set cj) = p (see later), we do not need to write down explicitly 
the constraints for this theory. One may easily recover the classical equations, i.e., without 
the trace anomaly term, by setting k = in the equations of the original theory. 

We consider first solutions of the classical equations. Equations (2.8b,c) have the 
solution 

e -2<i> = e -2 P= M_ x 2 x+x - j (212) 

A 
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where M is a constant, and arbitrary additive constants to x ± have been ignored. If M = 
we obtain the so-called linear dilaton, but if M 7^ the solution represents a black hole of 
mass M, with horizons given by x + x~ = and a singularity when x + x~ = M/A 3 . 

Consider next the situation where a linear dilaton occurs for x + < x£. At x^ a matter 
wave described by / = f(x + ), which is a solution of (2.8a) propagates in the ^"-direction. 
If f(x + ) has compact support, then once the wave has passed, the spacetime will once 
again be described by (2.12), but now we must expect M^O. For simplicity we consider 
an impulsive wave described by 

kJ2( d +f^ = a5 ( x+ - x o)^ ( 2 - 13 ) 

i 

where a is a constant. 

We now apply the constraint (2.11a) on an initial surface x~ = x~ . On such a surface 
the value of p is arbitrary; changes correspond to a rescaling of coordinates. We may 
therefore choose p = (f) on this surface. Then (2.11a), (2.13) imply a jump increase in 
d + {e- 24> ) at x + = x+, i.e., 

e -24> = e -2 P = ax + _ A 2 ^- + _^ x+ (214) 

for x + > x+ on x~ = x~ . Comparing this data with (2.12) we see that we have a black 
hole solution 

e -20 = e -2 P = ax + _ A 2 ^- + ^ (215) 

for x+ > x+. An alternative, more long-winded approach (in this case) is to solve (2.8b,c) 
as a characteristic initial value problem. If p and 4> are specified for x + > x+ on x~ = x~ , 
and for x~ > x~ on x + = then the solution is determined locally and uniquely for 
x ± > x ±0 . One piece of the data is given by (2.14). The other follows from continuity at 
x + = x+, viz. 

e~ 2 ^ = e ~ 2p = -\ 2 xtx~, (2.16) 

for x~ > x~° on x + = x„. 

We turn now to the semi classical analogue. For x + < x+ the classical solution is 

e- 2( t> = e ~ 2p = -X 2 x + x~. (2.17) 

This is also a solution of the semi classical equations. (Both sides of (2.8b) vanish for 
arbitrary k.) After the shock we have to apply the constraint (2.11a) for x + > x^ on 
x~ = x~ . We still have the coordinate freedom to choose p = cj) on this surface. We 
wish to study the situation in which there is no incoming energy momentum apart from 
the matter wave. This corresponds to choosing t + (x + ) so that the factor multiplying k in 
(2.11a) vanishes on x~ = x~ . Similarly, we want no energy momentum coming from the 
linear dilaton region. This corresponds to choosing t-(x~) so that the term multiplying k 
in (2.11b) is zero. With this choice equations (2.14) and (2.16) form characteristic initial 
data for the semi classical evolution equations (2.8b,c). However we lack an exact solution 
to these equations in the absence of some special symmetry. We shall therefore resort to 
numerical integration in section 4. 
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3. The conformal and conserved current theories 

In the conformal theory of de Alwis and Bilal & Callan the Lagrangian may be written as 



1 

L =- 

7T 



e~ Z(p (2d + pd-(f) + 2d-pd + (p - 4d + (pd-(p) + 
1 N 



i=i 



(3.1) 



where D(<f>) was given earlier by equation (2.6). Bilal and Callan suggested a sequence of 
changes of dependent variables 



uj = 



where e = k/\k\, followed by 



O = ^euj\/uj 2 — e — \ log(o; + a/cj 2 — e). 
These produce a free field Lagrangian 



7T 



and constraints 



1 N 

4Kd+n d-n - 4 K d +x d- X + -z d +h d -h + ^^ P ~ 24> D{4>) 



1=1 



(3.2a) 



(3.26) 



(3.3) 



1 N 

2nd ± 2 X + 4nd ± n d ± n - 4kO ±X d ±X + -J2 d ±& + *± ( a± ) = °- ( 3 - 4 ) 

1 i=i 

Note that Bilal and Callan used rescaled asymptotically Minkowskian coordinates 

a + = log(x+), a~ = - log(-aT), 

where x ± are the coordinates used in section 2. 

The equations of motion simplify after a further change of dependent variables ^± = 
X±0 to 

<9+d_tf_ = 0, d+8^ + = -^e 4 *-. (3.5) 

The first equation is the standard wave equation in characteristic coordinates, and the 
second is similar but with a known source term. Bilal and Callan wrote the solution in 
the form 



2^ + =2 1 (a + ) + 25(a-)-a(a + )-P(a-)+K~^J e 2a(s) ds J e 2m dt, 



(3.6) 



K = I + lo, 



2 \ A a/Ik 



where 



and a, (3, 7 and 5 are arbitrary functions of one variable, to be determined from the initial 
data and constraints. Bilal and Callan chose to take t±(a ± ) = in (3.4), which can be 
rewritten as 



(a_5) 2 - (d.(6 - (3)f - dJ6 =i- f>_/,) 2 . 



1=1 (3.7) 

^ AT v ; 



i=l 



The linear dilaton is not a solution of this theory. Consider however static solutions, 
i.e., depending on a = \(o + — a~) only, which are asymptotic to the linear dilaton as 
u — ► 00. Bilal and Callan obtained the general solution in the form 



a(a+) = ±a+, 7 (^ + ) = K + \T + \ log 



\k\ 

Te (3.8) 
P(a~) = -\a~, 8{<t-)=-\<t- + \T, 



where T is a constant that behaves like the mass. 

Bilal and Callan modelled the shockwave problem by requiring the solution (3.8) to 
hold for a+ < and setting \ E09_.fi) 2 = 0, \ Y.( d +fif = a5(a+) in (3.7). The solution 
is 

«(^ + ) = k + , 7(^ + ) =K + \T - V + - l)9(a+) + \ log M, (3 _ g) 
/3 ( (7 -) = _i (T -, %-) = - i<r- + |T, 
and (3.6) now implies that for a + > 

= - ^ + - 1) - - O + T + i log M, 

p + iogA + ie- 2 *^"*-"" -i)-i(<r+-0 + r. 

K K K 4 



(3.9) 



The variables % and O are regular functions of position. However there may be a 
singularity where the curvature scalar R = 8e~ 2p d + d-p diverges. In order to locate this 
recall from (3.2) that p = 1\ — ecu 2 and O = Q(uj). Then 

d+d-p = 2d+d- X - - ( 1 - „>^ J d+ft d-Q. 

Thus the singularity occurs when O' = 0. However from (3.2b) we see that this occurs 
when uj 2 = e or O = and only for k > 0. The apparent horizon is located where d+<fi = 
or equivalent ly where d+Vl = 0. 
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We now demonstrate that in this theory the singularity is eventually naked, i.e., the 
apparent horizon moves to the future of the singularity. Let the singularity be located at 
a~ = a~(a + ), and the apparent horizon at a~ = a^(a + ). If T is positive, the singularity 
will start off inside the apparent horizon. However for large a + (3.9) implies 

cr~ ~ — log(a + 7Kff + e~ ff+ ), 

SV 4 _ + ' (3.10) 

a h ~ - log(a+ \ne~ a ). 

In this limit a~ approaches from below, i.e., the apparent horizon and the singularity 
are ultimately tangent with the singularity to the past of the horizon. So the apparent 
horizon and singularity must meet and cross at some point (cr+, cr+). See figure 1. At 
this point the singularity will become timelike and naked. The line a~ = a~ , a + > cr+ 
will become a Cauchy horizon. Although the exact solution continues smoothly beyond 
the Cauchy horizon, it is unphysical because it has a steady outflow of radiation and an 
effective mass parameter that becomes arbitrarily negative. It should be noted that the 
argument does not depend on the impulsive nature of the shock, and can be generalized 
easily to arbitrary but finite infalls of matter. 

The analysis of the conserved current theory of Russo, Susskind and Thorlacius is 
very similar. They use the coordinates x^ of section 2 and auxiliary variables 

n = ^+ x = v^(p-0) + n. (3.ii) 



The Lagrangian is 
1 



S = 

TV 



(3.12) 



d + n a_o - d +x d- X + A 2 e 2 ^-")/^ - 1 d+fi d.u 

L i=l 

with constraints 

1 N 

V^d 2 ± x - d ±X d ±X + d ± n d ± n + -J2 9±fi d±fi - Kt±(x±) = 0. (3.13) 



i=l 



The asymptotically flat static geometries with <j) = p are given by 

n = X = - + Pv^log(-A 2 x+x-) + ^L, (3.14) 



where P and M are constants. Setting P = — j and M — gives the linear dilaton vacuum. 
Russo et al constructed a solution which is the dilaton for x + < x+ and corresponds to 
infalling matter for x + > x^ with 



1 N 

-J2(d+fi) 2 = a5(x+-xt), 



i=l 
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VIZ. 

n = X = - - ^ log(-A 2 x+x") --£=(*+- *+)0(* + - *+). (3.15) 

The singularity is again given (for k > 0) by O' = which implies O = ^y^l — log 
Again the singularity starts off inside the apparent horizon. However for large x + , the 
singularity will be at 

a k logx + 
Xs ~~ A 2 ~ lA 2 ^^' 
while the apparent horizon is located at 

a k 1 



h ~~ A 2 4A 2 x+ ' 

Again the singularity and apparent horizon must meet and cross at a point (cc+, x~) where 
the singularity will become timelike and naked. In this case one can cut off the solution 
(3.15) on the Cauchy horizon x~ = x~ and join on to the linear dilaton solution. This 
makes the Cauchy horizon a mild thunderbolt singularity. 



4. Numerical results for theories without special symmetries 

In both the classical and semi classical problems a numerical treatment is not entirely 
straightforward, for singularities are present in each. Fortunately we know the analytic 
solution for the former. We have therefore developed a numerical algorithm which handles 
the classical problem in a satisfactory manner, and have then applied it to the semi classical 
case. 

The parameter A may be scaled away, and so we have chosen A = 1. The values 
defining the initial data surfaces are arbitrary, and we have chosen x^ = ±1, except 
where stated otherwise. From (2.15) we see that the black hole singularity occurs on the 
hyperbola 

*--£(£-!)• («) 

and the relevant apparent horizon is at 

x- = (4.2) 

Thus by choosing a G (0, 1) we may ensure that the singularity remains in the domain of 
dependence of our data. The graphs are all drawn for a = 0.9, thought to be a generic 
case. 

At the singularity the Ricci curvature scalar R = 8e~ 2p d+d-p becomes singular. 
Figure 2 shows arctani? as a function of x ± for the classical case. The solution is not 
defined to the future of the singularity (top right of the surface) but for convenience in 
drawing the surface R has been assigned a token value of oo in this region. Also shown 
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is the apparent horizon x = —a. For large x + the singularity approaches the apparent 
horizon as predicted by (4.1,2). 

The same computational algorithm was adopted for the semi classical equations. From 
(2.8b) R may be expected to become singular once <p has increased to a value <p c at which 
P = 0. Such a 4> c exists for the original theory if k > and for the decoupled ghost theory 
if R > 2. We have also considered k < 2 that is iV < 24 for the decoupled ghost theory. For 
each fixed x + the programme integrated the equations in the direction of increasing x~ . 
The singularity was deemed to occur at the point where = <p c and the apparent horizon 
when d+cj) changed sign. For k > 0, the results do not seem to depend sensitively on the 
exact value so k = 0.5 or 0.8 was used for the graphs presented here. The behaviour shown 
in Figure 3 for the original theory looked at first sight broadly similar to the black hole 
case, although the singularity is a little steeper. However the two solutions are radically 
different. To show this we need a test that will distinguish a black hole singularity that 
remains at a fixed position from a thunderbolt that spreads out to infinity. 

Consider the outgoing null geodesic x~ = x~ , with tangent vector T + = dx + / dt where 
t is an affine parameter. For each fixed value of x + we consider the ingoing null geodesic 
x + = const, with tangent vector T~ = dx~ / ds and affine parameter s normalized by say 
g(T + ,T~) = | at (x + ,x~). The affine parameter distance to the horizon will be denoted 
s(x + ). For a black hole we would expect that for large x + , s(x + ) will be asymptotically 
a linear function of t(x + ), c.f., Schwarzschild. Indeed it is an elementary exercise to carry 
out the calculations analytically for the black hole in the classical shockwave problem, 
finding 

t( x +) = - (a + x-)- 1 log (ax+ - (a + x~)x+) , 

s(x + ) =(-ax+/x + + a + x~) [log(ax+) + (a + x+)t(x + )] . 

As x + — > oo, t — > oo and s(x + ) ~ (a + x~) 2 t, giving linear behaviour. 

We next developed a numerical algorithm to explore the behaviour of s(t) for large t. 
Numerically this is not entirely straightforward. Firstly one wants to explore very large 
values of x + , i.e., to integrate over an enormous number of grid points. Secondly the 
singularity is approaching the horizon asymptotically, c.f., (4.1,2)! Figure 4 shows the the 
computed and analytic behaviour of s as a function of t for the black hole arising in the 
classical shockwave problem, demonstrating the stability and accuracy of the algorithm. 
As expected the behaviour is asymptotically linear. 

When the same algorithm is applied to the solution of the semi classical equations 
significantly different behaviour is encountered. As can be seen in Figure 5, s(t) is definitely 
non linear and appears to be either bounded above, or at most logarithmic. This behaviour 
is also observed for other values of the parameters. We take this as an indication that the 
singularity does not remain in a bounded region, like in a classical black hole, but spreads 
out to infinity thunderbolt. 

The decoupled ghost theory is also difficult to treat analytically and so we resorted to 
numerical computation for this theory as well. One of the arguments that were advanced 
for this theory was that coupling the ghosts to a different metric would mean that a black 
hole wouldn't radiate a negative energy flux if the number N of scalars was less than 24. 
We therefore tried iV = 12, k = —1. However the numerical results shown in figures 6 and 
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7 are radically different. We interpret them as appearing to indicate a black hole that is 
growing in size with the apparent horizon moving out. Presumably this implies that the 
energy flux of the outgoing radiation is negative. We expect this to be true in any of the 
four theories if k is negative. We therefore calculated the more physically reasonable case 
with k = +1. This was similar to the original CGHS theory. Figure 8 shows arctani? as a 
function of x ± . As in the earlier cases there is a singularity which is located asymptotically 
at x~ = const, as x + — > oo. The apparent horizon lies to the past of the singularity and 
appears tobe asymptotic to it. The behaviour of t(s) shown in figure 9 is also similar. We 
had to integrate a lot further in this case but again it appeared to be bounded indicating 
a thunderbolt singularity. 

5. Conclusions 

If an evaporating black hole does not settle down to a stable remnant, any attempt to 
describe it by a Lorentz metric must have either a naked singularity or a thunderbolt. We 
have studied four toy two dimensional models of black hole formation and evaporation. In 
the two theories whose Lagrangians had extra symmetry, the conformally invariant and 
conserved current theories, it was possible to write the general solution in terms of new 
variables, x an d ^- The solution was non singular in terms of these variables, but it had a 
naked singularity in terms of the physical variables, 4> and p. One might expect the semi 
classical approximation to break down at this singularity and not to determine the fields 
beyond the Cauchy horizon that starts at the point where the singularity first becomes 
visible from infinity. 

In the case of the two Lagrangians with extra symmetry, the Cauchy horizon was 
regular when approached from below. However we suspect that this will not be the case for 
more general Lagrangians: we expect the singularity will be a spacelike or null thunderbolt 
that spreads out to infinity and means that spacetime can be evolved only a finite retarded 
time according to the semi classical equations. This expectation is strengthened by the 
numerical calculations we have done for two Lagrangians without special symmetries, the 
original model proposed by CGHS, and the decoupled ghost modification proposed by 
Strominger. The results point to thunderbolts in both cases if k is positive, i.e. if the 
number N of minimal scalar fields is greater than 24. 

One can interpret these results as follows. In two dimensions the conservation and 
trace anomaly equations seem to imply that a solution asymptotic to the linear dilaton 
will continue to radiate at a steady rate. Either the evolution of the solution will be cut off 
after a finite time by a thunderbolt or the Bondi mass will become negative eventually. In 
the latter case, one might expect the singularity to change from being spacelike to timelike 
and naked on about the outgoing null line on which the Bondi mass becomes negative. The 
four theories we have considered illustrate the two possibilities: the two with additional 
symmetry give naked singularities with the Bondi mass becoming arbitrarily negative while 
the two more general theories give thunderbolts. 

We would expect a semi classical treatment in four dimensions (if it is possible to incor- 
porate back reaction consistently) would be similar to the more general two dimensional 
theories and would also lead to thunderbolts in general. In our opinion, the prediction 
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of a thunderbolt would indicate not that spacetime came to an end when a black hole 
evaporated, but that the semi classical approximation broke down at the end point. We 
would expect a full quantum treatment would soften the thunderbolt singularity into a 
burst of high energy particles. It would be tempting to try to connect such events with 
the gamma ray bursts that have been observed, but there is a problem with the energies 
involved. There is no reason to expect the semi classical approximation to break down 
until the horizon size becomes of the order the effective Planck length. In the case of black 
holes without a conserved charge, this will not happen until the black hole gets down to 
the Planck mass, so there's far too little energy left to explain the observed gamma ray 
bursts, specially if they are at cosmological distances, as the observations seem to indicate. 
Black holes with a conserved charge but in theories without a dilaton field approach a zero 
temperature extreme state, so the semi classical approximation shouldn't break down and 
there's no reason to expect a thunderbolt. In theory with a dilaton field with the coupling 
to gauge fields suggested by string theory, the semi classical approximation can break down 
while the black hole still has a macroscopic mass. However, the mass difference between 
the black hole at this point and the zero temperature extreme black hole, which is presum- 
ably the ground state with the given charge, is much less than the Planck mass. In fact it 
is even less than one quantum at the temperature of the black hole. So any thunderbolt 
predicted by the semi classical approximation would have to be extremely mild and could 
not account for the observed gamma ray bursts. If the universe does contain black holes 
that are reaching the end points of their evaporation, it seems they will do it without much 
display. 
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Appendix. Numerical Methods 

This would seem to be the first paper in this area to utilize explicit numerical solutions 
of the field equations. This is somewhat surprising for in two dimensions reliable accurate 
numerical solutions can be obtained readily using even modest computer workstations. 
The purpose of this section is to explain in some detail how our numerical solutions were 
obtained, so that our methods become accessible to others. 

The fields are functions on the x^ plane. We replace the plane by a two dimensional 
lattice with equal spacing h in the x + and x~ directions. Figure 10 shows a typical lattice 
cell. The four corners are denoted n, e, s and w, while the centre is denoted o. Let 
y(x + ,x~) be a function taking values in R n and let y n , y e , y s , y w and y Q be the values at 
the corresponding grid points. We assume that the function y is sufficiently regular that 
it can be represented within the cell by a Taylor series with remainder term 0(h 4 ). It is 
then a routine exercise to verify the following relations: 

y =^(y w + y e ) + 0(h 2 ), (Ala) 
yo=l(yn + y e + y s + y w ) + 0(h 2 ), (Alb) 

(d + y) o =^^+0(h), (A2a) 
(d + y) = V°-V'+V»-V™ + (h 2 ), (A2b) 
together with the obvious analogues for (<9_y) D , and 

(d + d_y) o = (y"-ye + y*-y*>) + ( h 2 ). (as) 

All of the theories treated here have field equations of the form 

d+d-y = F(y,d ± y) (A4) 

where y = (p, 0) and F is smooth. Further, initial data is given on the initial surfaces 
x^ = xf. If we discretize (A4) and the initial data according to the above prescription 
then the paradigm problem is the following: given y s , y e and y w , determine y n . Evaluating 
(A4) at the point o and using the relation (A3) we obtain 

Vn =y w + y e -y s + h 2 (F( y , d±y)) + 0(h 4 ) 

=y w + y e -y s + h 2 F ( Vo , (0 ±y ) o ) + 0(h 4 ), (A5) 

where the last transformation is a tautology. Nevertheless (A5) is the basis for our numeri- 
cal algorithm. We propose to evaluate it iteratively twice, leaving y w , y e and y s unaltered, 
but replacing the arguments of the function F by approximate values. 

For our first evaluation we use approximations (Ala), (A2a) finding 

y n := y e + y w - y s + h 2 F (\(y w + y e ), hr x (y e - y s )) + 0(h 3 ). (A6) 
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Since all of the explicit terms on the right hand side of the equation are known we may 
evaluate a trial approximation to y n , and we have used here an atom of PASCAL formalism, 
":=" , whereby the evaluated right hand side of the equation is then assigned to the labelled 
quantity on the left hand side. 

We can however do significantly better than this. For our second evaluation we use 
approximations (Alb), (A2b) finding 

Vn ■= y e + y w -y s + h 2 F (\{y a + y e + y s + y w ), ^h~ 1 (y n -y w + y e - y s )) + 0(h 4 ). (A7) 

In principle we could repeat this, regarding it as an iterative process for solving the non- 
linear equation (A7) for y n . However subsequent corrections to y n are smaller than the 
truncation error (9(/i 4 ) inherent in the equation. Although the improvement in the error 
bound of (A7) over (A6) may look small it is essential. In order to integrate the equations 
out to large x + the algorithm has to be applied ~ h~ 3 times, and the errors committed at 
each stage are cumulative. 
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Figure Captions 



Figure 1. 

This figure shows some features of the Bilal and Callan exact solution for a = 1.0, k = 1.0 
and T = 4.0. The four curves show the positions in the x^ plane of the singularity, the 
apparent horizon, what Bilal and Callan call a "horizon" and the Cauchy horizon. 

Figure 2. 

The surface drawn is z = arctani?(x + , x~) where R is the Ricci curvature for a classical 
black hole with a = 0.9. The solution is not defined to the future (right) of the singularity, 
and so R is assigned a token value of oo. Also shown is the apparent singularity, where 
d + 4> changes sign. 

Figure 3. 

The surface drawn is z = arctani?(x + , x~) where R is the Ricci curvature for the original 
CGHS theory with a = 0.9, k = 0.5, N = 30. The solution is not defined to the future 
(right) of the singularity, and so R is assigned a token value of oo. Also shown is the 
apparent singularity, where d + (j) changes sign. 

Figure 4. 

The affine parameter distance s along an ingoing null geodesic from the initial surface 
x~ = x~ to the apparent horizon is plotted against t, the affine parameter distance along 
the initial surface. Both curves refer to a classical black hole with a — 0.9. The solid line 
was obtained by numerical integration, while the dashed line was computed analytically 
from equations (4.3). As t — > oo the relation becomes linear. 

Figure 5. 

The affine parameter distance s along an ingoing null geodesic from the initial surface 
x~ = x~ to the apparent horizon is plotted against t, the affine parameter distance along 
the initial surface. The solid line was obtained by numerical integration of the original 
CGHS theory, while the dashed line was computed for a classical black hole with the same 
initial data. In both cases a = 0.9, and for the CGHS theory k = 0.8, A = 30. As t — > oo 
s appears to be bounded above. 

Figure 6. 

The surface drawn is z = arctan R(x + , x~) where R is the Ricci curvature for the decoupled 
ghosts theory with a = 0.9, k = —1, N = 12. The solution is not defined to the future 
(right) of the singularity, and so R is assigned a token value of oo. Also shown is the 
apparent singularity, where d + (j) changes sign. 
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Figure 7. 

The affine parameter distance s along an ingoing null geodesic from the initial surface 
x~ = x~ to the apparent horizon is plotted against t, the affine parameter distance along 
the initial surface. The solid line was obtained by numerical integration of the decoupled 
ghosts theory, while the dashed line was computed for a classical black hole with the same 
initial data. In both cases a = 0.9, and for the semi-classical theory k = — 1, N = 12. As 
t — > oo s appears to be unbounded above. 

Figure 8. 

The surface drawn is z = arctan R(x + , x~) where R is the Ricci curvature for the decoupled 
ghosts theory with a = 0.9, k = 1, N = 36 and x+ = 4. The solution is not defined to the 
future (right) of the singularity, and so R is assigned a token value of oo. Also shown is 
the apparent singularity, where d+4> changes sign. 

Figure 9. 

The affine parameter distance s along an ingoing null geodesic from the initial surface 
x~ = x~ to the apparent horizon is plotted against t, the affine parameter distance along 
the initial surface. The solid line was obtained by numerical integration of the decoupled 
ghosts theory, while the dashed line was computed for a classical black hole with the same 
initial data. In both cases a = 0.9, and for the semi-classical theory k = 1, N = 36, 
x+ = 4. As t — > oo s appears to be bounded above. 

Figure 10. 

The computational grid in the x ± -plane. The plane is replaced by a lattice with spacing h. 
Given data at points s, w and e, the numerical algorithm estimates d+d-y at the fictitious 
point o and hence the dependent variable y at the new lattice point n. 
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